# Clean up
rm(list = ls())

# Set working directory; please set your own here
setwd("~/Dropbox/Fragmentation/fragmentation_replication_bjps/")

library(tidyverse)
library(readxl)

# Import the election-level dataset
first_stage <- read.dta13("Data/fragmentation_electionlevel.dta",generate.factors=TRUE)

# Import data with main analyses
data <- read.dta("Results/analyses_final.dta")

# Creating a theme for the plot
theme_base <- 
  theme_minimal(base_size=9)  + 
  theme(legend.position = c(0.85, 0.1),
        axis.text=element_text(size=12),axis.title.x=element_text(size=12),axis.title.y=element_text(size=12),
        plot.title = element_text(size=12, hjust= 0.5),
        legend.text = element_text(size=12))

# Plotting the total number of parliamentary parties against the added number of parties by parties in the bandwidth
partiestotal_plot <- ggplot(first_stage, aes(x = inbw50, y = cnt_parl_parties)) +
  geom_jitter(width = 0.15, alpha = 0.6) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red",size = .5) +
  stat_smooth(method = "lm", aes(colour = "Linear model")) +
  stat_smooth(method = "loess", linetype = "longdash", aes(colour = "LOESS")) +
  scale_colour_manual(name = " ", values = c("blue", "green3")) +
  labs(x = "Total number of parties within the bandwidth", y = "Total number of parliamentary parties") +
  scale_y_continuous(limits = c(-3, 18)) +
  #  annotate("text", x = 2, y = -2, label = "F(1, 35) =  17.45", size = 4) +
  theme_base 
partiestotal_plot

# Plot distribution of f-statistics
firststage_distr_total <- data %>%
  filter(Sample == "Whole") %>%
  filter(FE == "None") %>%
  filter(id == "No. parties just above") %>%
  filter(Predictor == "Absolute number parl parties") %>%
  ggplot +
  geom_density(aes (x = fstat), fill = "steelblue", color = "blue", alpha = 0.5) +
  scale_x_continuous(name = "F-test", limits = c(0, 20)) +
  scale_y_continuous(name = "Density") +
  geom_vline(aes(xintercept = mean(fstat)) ,col = 'red',size = 1, linetype = "dashed") +
  theme_bw() 
firststage_distr_total

# Put the two plots together
first_stage_totalparties <- grid.arrange(partiestotal_plot, firststage_distr_total, ncol = 2, nrow = 1)

# Save final plot
ggsave("Plots/figa5.png", plot = first_stage_totalparties, width = 30, height = 20, units = "cm")